{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([[ 1.      ,  5.      , 21.944048],\n",
       "        [ 1.      ,  4.      , 22.171763],\n",
       "        [ 1.      ,  2.      , 22.440421],\n",
       "        [ 1.      ,  4.      , 22.457996],\n",
       "        [ 1.      ,  4.      , 27.728212]]),\n",
       " (400, 3))"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "\n",
    "#加载数据\n",
    "def load_data(file_name):\n",
    "    with open(file_name) as fr:\n",
    "        lines = fr.readlines()\n",
    "\n",
    "    data = np.empty((len(lines), 3), dtype=float)\n",
    "\n",
    "    for i in range(len(lines)):\n",
    "        line = lines[i].split('\\t')\n",
    "        data[i] = line\n",
    "\n",
    "    #根据y排序,主要是为了画图\n",
    "    argsort = data.argsort(axis=0)[:, 2]\n",
    "    data = data[argsort]\n",
    "\n",
    "    return data\n",
    "\n",
    "\n",
    "data = load_data('自行车速.txt')\n",
    "data[:5], data.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#定义常量\n",
    "N = len(data)\n",
    "\n",
    "M = data.shape[1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3dcYxc1ZXn8e/pdjtjZ7IxXpwMNDYGi3EGx8GOW8DKu1miTMYkQtCBZYIniRhttE4kop1oIiswaylmRRS0iGT/mJ2MiILIisQhwdABkV0SkYzIRDHQxhjjAJvA2MaNgz0YQxQbMO2zf1SVqWq/97rfedWvql79PlKru291Uc/Vzalb5557rrk7IiJSLQOdvgAREWk/BXcRkQpScBcRqSAFdxGRClJwFxGpoDmdvgCA008/3ZcuXdrpyxAR6Snbt2//V3dflHRbVwT3pUuXMj4+3unLEBHpKWa2N+02pWVERCpIwV1EpIIU3EVEKkjBXUSkghTcRUQqqCuqZURE+s3YjgluefBZXjxyjDMXzGPjuuWMrh5u239fwV1EpGRjOya44Z5dHDs+CcDEkWPccM8ugLYFeKVlRERKdsuDz54M7A3Hjk9yy4PPtu0xFNxFREo2ceRYrvEIBXcRkZINmuUaj1DOXUSkoLyLo5MpJ+CljUcouItIJc12NUrz4+RdHB1eMC8xBTO8YF7brktpGRGpnEbAnThyDOftgDu2Y6LtjxVZHN24bjlDg60pmKFBY+O65W27LgV3EamcMqpRGsKLo1MzMO3LyAAK7iJSQS+mBNa08SIii6O3PPgsx0+0RvPjJ1ylkCIiWc5MyV2njRcRWRxVKaSISMDGdcuZNzTYMjZvaLCtOe2GyMx9IOWmtPEIVcuISOU0qlTKqJaJzNxPpNyUNh6h4C4ilTS6enhWgvlUg2aJgbydG5IilJYRESkgMnNPC/vtfDnQzF1EKqmsTUwDlpxOycqfp4X9dlZDauYuIpUztmOCjXfvbNnEtPHunbOyiSmSP0/biVrqDlUzu93MDprZU01jd5nZE/WPPWb2RH18qZkda7rtH9t2pSIiM3Tj/bs5PjmljnzSufH+3R26olYfft+iXOMRM0nL3AH8PfC/GwPu/snG12Z2K/Bq088/5+6r2nWBIiJ5vXL0eK7xIhbMG+LIsVP/uwvmDaXe5+fPHMo1HjHtzN3dHwYOJ91mZgb8JbClbVckItJDNl++4pT8+oDVxtOUsYO26ILqfwBecvffNI2dY2Y7gNeATe7+i6Q7mtkGYAPAkiVLCl6GiEhx0UXYQTNONFXHTFcGeWZKV8h27qAtuqC6ntZZ+wFgibuvBv4W+J6Z/ZukO7r7be4+4u4jixa1L88kItUztmOCtTf/jHOuf4C1N/9sVhZGo50kI31iysi5h4O7mc0BrgTuaoy5+xvu/nL96+3Ac8CfFr1IEelfZbXvjXaSjKRYuiLnnuHPgWfcfX9jwMwWmdlg/etzgfOA54tdooj0s7La90bz4JEmZWXk3GdSCrkF+BWw3Mz2m9ln6zddw6kLqR8CnjSzncDdwOfdPXExVkRkJspq3xvtJBlJsSyYn1xJkzYeMe2CqruvTxn/64SxrcDW4pclIlJTxuIj1ILxndv2JY5niaRY0joTtPEIVe1QFZFy5V0cLWPxEeCBJw/kGm+IvLN4NaEuPms8QsFdREoTaQsQmRlHeqxHNz5F0jllpGUU3EWkNJG2AJGZ8cXnnpZrvIjIwSBKy4hIpURmx5FZ7p6XkwN/2ngRo6uH+eCSd7eMfXDJuzM3PyktIyJ9740pZZDTjUO5B2RvGtvFL59rLQr85XOH2TS2K/U+ZZzxquAuImFl7Bw9evxErnEoJ6fdsOWRF3KNQ5fvUBWR/hbZOZrWKTGrg2JEZLYfPbQ6chJTtDInDwV3EQmJ7BzdfPmKU46SM7I7KEaOpIvM9t8xJzkcpo0XUUZLYgV3EQmJ5LXH9x4+5Sg5r4+nKeNIOoDXUwJ/2ni3U3AXkZDIomAkPx2pWY/M9svM05eRnlJwF5HQwmikvjuSn47cJzLbj+Tpoy674Ixc4xEK7iJ9LtpSd3T1MFetGT45gx4046o1w5n13ZFZeOQw6dNSZttp4xDL00NsFt7tLX9FpAKiLXXHdkxw12MvnJxBT7pz12MvZL4orL9oca5xqL1DSDrGLusdQpmz8G49Zk/BXaTPRQNNpJXATaMr+fTFS1pm+5++eAk3ja5Mvc8Px/cx5aAjTnhtPE2oNr5AHnzqO4+ZHLOXZzxCwV2kz0UXEqPlfCNnL+RP3v1HGPAn7/4jRs5emPnzU3d/TjcetfnyFQxNmYIPDVjmDBxix+xF1ivyKnpAtoh0mbyHPJfRxKr52jbevfPkjL/RFRKY0UHUM2UkL55mzacbj5/3gOykXvNZ443HGt97mC2P1NJaM1mvyEvBXaRCGoujjRx6Y3EU0oNntInVgnlDHEn4maw0RlYqp52BLVobP7o6f4AdNEus3MlKzYztmGDr9omW9Yqt2ycYOXth254HpWVEKiSyOBrN/0Z2m5axMxNiFTZRkVLNMs6FVXAXqZDI4mg0/xvZbRoRKWss6/QmiL2QqFpGRHKJzMIj9eoQ220aEVkTKKMxV0PkxbErqmXM7HYzO2hmTzWNbTazCTN7ov7x8abbbjCz35rZs2a2rm1XKiLT2rhuOUODUyo+Bi0z0ETq1Rs/l2c8KrImUFb6B2Ivjt3S8vcO4NKE8W+4+6r6x48BzOx84BpgRf0+/2Bmgwn3FZHZkpQryRCpV4fYbtNILfn8uckhJG28bGmLo+0+FzavaYO7uz8MzDSJdgXwfXd/w93/BfgtcGGB6xORHCI119FZbmS3aaSW/OibybtK08bLFlkc7fac+xfM7Ml62qZx6uww0Jxw218fO4WZbTCzcTMbP3Sofa9WIv0sUnMdFdltOrp6mE9euLjlPp+8cHF2HX7OcSjvUBCIBequyLmn+CawDFgFHABurY8nvR9L/B24+23uPuLuI4sWtX8FW6QfRU4TKhII8+42jaQwIumf6G7TSHfMSKAuY4dqKLi7+0vuPunuJ4Bv8XbqZT/Q/J7sLODFYpco0r82je1i2Q0/Zun1D7Dshh9nHroMnNKDZbpxKBYI83aTjKQwzl00P9c41N4h3HL1BQwvmIdRK0u85eoLMt8hRLtjRgL16Ophvnblypbr+9qVKzu/Q9XMznD3Rk3RJ4BGJc19wPfM7OvAmcB5wKOFr1KkD20a28Wd295ujjXpfvL7rNRHXklb4adLlUB2oE67bySF8fyho7nGG/LuNo38exqP07h/nrYFkd2weUwb3M1sC3AJcLqZ7Qe+AlxiZquopVz2AJ8DcPfdZvYD4NfAW8B17t4dqx4iPSarjrydwT2tFHK6rfDRXHNS/j8rhREtudw0tqvlBWv9RYszn7cii5yzHagjZlIts97dz3D3IXc/y92/7e6fcfeV7v4Bd7+8aRaPu3/V3Ze5+3J3/z+ze/ki1RUJavOHkv+XThuHeClkpJtkpA4/knNvvOtpfsG6c9u+zLRWGYucZdIOVZGS5M2fR8ydk1z7nTYO8VLIcDfJnHX4kZLLyO7ZMhY5y6TgLlKCyEwyItrhsazHitThj5y98JRANVAfTxN511PGImeZFNxFShCZSUbSEZHUQrQUMpKWieS1b3nwWaaen3SiPp4m8txVjYK7SAkiM8nIfSI57RVnvivXeEMkLRN58Ym8IERSOdFSyG6l4C5SgshM8p0pvVPSxhsmp6Q9pn4/1bbnX8k13hBJy0Ty2pF3CJHds2X0WC+TTmISKcH6ixa31Kw3j6f5Q0rvlLRxqFW+JB0mnXXSUbTUcMH8ocRF16ygO7p6mB+O72s5//SDS949K8cA3jS6MlfJaBn9XsqkmbtICW4aXcnaZa0LgGuXLWxrvTqU2+o2EnQ3je065WDrXz53OHNhuaxFYpVCikhuYzsmeHzfqy1jj+97tWfzuUDi+alZ4xBbWI6kZSJUCikiuZWVz507mJzDTxsvItKkLJICeuN4choqbTyqaqWQyrmLlKCsfO6bk8lBMm28iEiTsgFLvj3rBeHo8amFkNnjRXRjG4EozdxFSlC1fG7UO+Ykh5y0cYnTMypSgsiZmfNS+sGkjfeC11Nm22njUO7BG1XSu38lIj0kcmbmVWvOyjUeFd3NmXZr1r0i72Ci/eb7nYK7SE6R03oix9+VcYgyxOvcP3XxklzjED/YIu/BG6IFVZFcGlvUG5UvjS3qQGawGTRLDJZZs+OyFmEj1wZvHxiSp2d6tx5sUUUK7iI5RE/ricyOIwdbRERn7pB/FygoUJdFaRmRHCLpFailEvKMQy2FMTgl1zw4kN0ELCJybdL9FNylr+XNn0cXHyO55vG9hxObgI3vPZxyj5iq7cyUGgV36VuRFq/RFEZk92Nkq37E6Ophrloz3NJB8ao1s5c6iSxIS37KuUvfiuTPh1Py4DNJYeTNNUdeSE5L6dR4WkYflrEdE2zdPtFyStTW7RPTHpAdEV2Qlvw0c5e+FalGiRyGERVJAUU6NZbZx7xqPdO72bTB3cxuN7ODZvZU09gtZvaMmT1pZvea2YL6+FIzO2ZmT9Q//nE2L16kiHBLgJwHPDfkTUdEThOKdGoss4951Xqmd7OZzNzvAC6dMvZT4P3u/gHg/wE3NN32nLuvqn98vj2XKdJ+kZYAkQOeIZbfj5wmFFFm3xv12CnPtDl3d3/YzJZOGftJ07fbgP/U3ssSmX2RHaDRUshofXykjjyvjeuWt+TBYebVMmM7JnJtSCryWJJPOxZU/zNwV9P355jZDuA1YJO7/yLpTma2AdgAsGRJ+nZlkdkSSRFEWtZGHysicn3RXaORxdHoY0l+hYK7mf034C3gu/WhA8ASd3/ZzNYAY2a2wt1fm3pfd78NuA1gZGSk/c2mRaYROQM00sMcyttt+o45AxxL6LA4XUvdyK7R6LsR7VAtR7haxsyuBS4DPuVeW4t39zfc/eX619uB54A/bceFirRb9ODliGiVTd5F2KTAnjVehBZHu1to5m5mlwJfBv6jux9tGl8EHHb3STM7FzgPeL4tVyoyjbz537IOXm5I2m2aJZL2iDYBiyjr3YjEzKQUcgvwK2C5me03s88Cfw+8C/jplJLHDwFPmtlO4G7g8+7e3r3SIgki1ShlHbwMcOP9u09J3Zzw2niaSE14kSZgealtQXebSbXM+oThb6f87FZga9GLEskrkv8tMy2TlNvPGodY2qPIDtq8tDja3dR+QCohEgjLTsvkFUl7lF1qqMXR7qX2A1IJkc0xZW6omZ9y7mnaOMRPLcrboEyqSTN3qYTIjPXD71vEndv2JY6329w5gxxNqFiZO2cw4adrdGqRFKHgLpUwunqY8b2HW458m65t7b2PJy+23vv4ROqu0PPe805+c/APieNZIj1fQIFa4hTcpSvlLWuMtK39w5uTucYBXjzyeq7xhjJLFEVAOXfpQpGyxrJayUZeEKDcEkURUHCXLhQJ1N2+W1LnlErZFNyl60QCdaTyZV5KpUraeBGR9sIiRSi4S9eJBOpI8PzalR845X+Agfp4u0XaC4sUoeAuXSdS3/3AkwdyjUOtEuWvphyG8VcXL8lcuF0wL6VlQcp4Q7enjaR6FNyl60Q24kS294/tmOCux15oqbC567EXMhduL7vgjFzjDTqBSMqmUkjpSmXUd994/26OT045Mm/SufH+3amPHU2v6AQiKZuCu1TC/KGBxB2gWdv7y2rmBWqyJeVTcJeulHcTU2R7f0SRHuaRdyN5nweRBuXcpeuM7Zhg4907WzYxbbx7Z2YuPNLhMbI4WmYP88hmLpEGBXfpOlm58DSRBcvI4ujo6mGuWjPcUmEzXQ+bqLJ23Uo1KbjLrMp7BijEcuGROvfI4mhaD5vZmE2rfFKKUHCXWVNmWiFS556UO88ah3Jn0yqflCIU3GXWlBkII7P9tI6MWZ0ay5xN64xSKULBXWZNt6cVIp0ay5xN61QlKWLaUkgzux24DDjo7u+vjy0E7gKWAnuAv3T3V+q33QB8FpgE/qu7PzgrVy5dr0jZYF4L5g0lHnyRVfkS6bGuM0qlV8xk5n4HcOmUseuBh9z9POCh+veY2fnANcCK+n3+wczaW2gsPaPMtMLmy1cwNNAalIcGjM2Xr0i9T2Tmrtm09IppZ+7u/rCZLZ0yfAVwSf3r7wD/BHy5Pv59d38D+Bcz+y1wIfCr9lyu9JIyd2VGHms45Z3FdD3WNZuWXhDdofpedz8A4O4HzOw99fFhYFvTz+2vj0mfKjMQju89zO9efR0Hfvfq64zvPZz52GUekC1StnYvqCYlKxPf45rZBjMbN7PxQ4fU01qK2TS2izu37WupP79z2z42je1KvY96rEuVRYP7S2Z2BkD988H6+H5gcdPPnQW8mPQfcPfb3H3E3UcWLdJMSYrZ8sgLucah+6t5RIqIBvf7gGvrX18L/Khp/Boze4eZnQOcBzxa7BKll0V2qEZ0e1mjSNlmUgq5hdri6elmth/4CnAz8AMz+yywD7gawN13m9kPgF8DbwHXuXv2sfDSM/J2KGzsUG2UDTZ2qAJtz8OXWdaoTo3SC2ZSLbM+5aaPpPz8V4GvFrko6T6RQJ21QzUrGEZ6s6+/aHHi4uj6ixYn/HRNpMKmzBcskSLUz11mJBKoy8xpj5y9kO9t20fzS8JAfTxL3mqe6AuWSNnUfkBmJNJka8H8lH7pKeMNSbP2rHGoBd2pt56oj7eTFmGlVyi4y4xEmmz94Y23co0XEXnxidAirPQKBXeZkUg1ypuTybeljRcxkPIakzYepU6N0isU3GVG0rbkT7dVvywnUl4v0saj1FtGekVPL6iqJK08ZXdD7GbqLSO9oGeDe+MQ5cZZm41DlEElabOhzCZgIlJczwb3rEOUFXBmR1kz1khvdhFp1bM598ixalJMWa0EIr3ZRaRVz87cpVyRnZlGckvQ6QpYlAISKU7BXWYksjMzrVBlJgUsWrQUKaZn0zJSrrI2CUWl5eOVp5d+1bPBPa2JVFZzKYmL7FAtk/L0Iq2UlpEZiexQHRqApHYws/H6qzy9SKueDe6R5lISFzlMOu1XMVu/IuXpRd6mHIbMSNqh0TpMWqQ79Wxw1wJauR548kCucUgveeyOLL1ItfVscL/sgjNyjUsxkU1jRUohRaSYng3uP3/mUK5xKV+3V9iIVFnPBvdur7uumkgaLFJhIyLt0bPBXbPCYvL2iYmkwfQ7EumccHA3s+Vm9kTTx2tm9kUz22xmE03jH2/nBTdoVhjX6BMzceQYztt9YrICfGRBVb8jkc4JB3d3f9bdV7n7KmANcBS4t37zNxq3ufuP23GhU2lWGJfVJyZNZEH1tJSDsNPGRaR92pWW+QjwnLvvbdN/b1qaFcaVtV6R9qvQr0hk9rUruF8DbGn6/gtm9qSZ3W5mpyXdwcw2mNm4mY0fOpS/wqWsA5El7tWEAzeyxkWkfQoHdzObC1wO/LA+9E1gGbAKOADcmnQ/d7/N3UfcfWTRovy7HMs6EFnizkxpTZA23qysg0FEqqodM/ePAY+7+0sA7v6Su0+6+wngW8CFbXgM6UEb1y1n3tBgy9hMDtWOLPiKSKt2BPf1NKVkzKy5Nu4TwFNteIxTaGt73GDKk5Q2HjW6epivXbmS4QXzMGpNxr525cppm3tFFnxFpFWhrpBmNh/4KPC5puH/YWarqO0y3zPltrbR1va4SEpr/tBAYsfN6frnRzo1vpiysJs2LiKnKhTc3f0o8G+njH2m0BXN0GnzhxLL8FRmN73IC+PcOYOJwX3unMGEny7mzJT2wjPJ1YtITc/uUFWZXbnKrHyJ5upF5G09e1hHNNiM7ZjQaT0BZc6mdaqSSHE9G9wjwaZRhdFYrGtUYQAKHNPYuG45X/rhTiabEvODAzZrs2mdqiRSTM+mZSJv3VWFUZO2Bpq1Njq+93BLYAeYPOGM7z3cxisTkXbp2eAeKbNTFUbN0GDyrz1tHGDLIy/kGheRzurZtAzkf+te1SqMvOsIkcPF1ctHpLf07Mw9oopVGGM7Jth4986W3Zwb796p3Zwifa6vgnt0x2Q3u/H+3RyfbJ09H590brx/d4euSES6QV8Fd6gtDP7u1ddx4Hevvt7zC4KRPusiUn09nXPPa9PYLu7ctu/k95PuJ7+/aXRlpy5LRKTt+mrmXsWKj8jB1d38OCLSHn0V3KtY8bHizHflGo/afPmKUw5CGbDauIh0n74K7lW07flXco0DrF22MNd4w9TzaXVerUj3UnDvcZF3I3teTt60lTYOtd29x6fsUD1+wvtud69Ir+ir4D6cslkpbbyqIgdka3evSG/pq+BexU1MEWnplKw0S5HzUEWkfH0V3EdXD3PVmuGTQWzQjKvWzE73wcgBz2UdCh1J5eiFUaS39FWd+9iOCbZunzgZxCbd2bp9gpGzF2YG+Ly9WyKthctsRzyc0mMnKz2lHusivaWvgntWy992Bt3I40TuA7U68yMJB5Rk1Z9/+H2LWjZzNY9nUY91kd7RV2mZyEJipAd8ZPExumC5+fIVTM2UG9n15z9/5lCucRHpPX0V3KduwpluHGJBN7L4GF2wHN97+JSDrb0+nib6QlLWmoCIFFcouJvZHjPbZWZPmNl4fWyhmf3UzH5T/3xaey61uBMp64Vp4xALumnpjay0R3TB8nuPnJpeyRqH2L+pkZ5qbi18wz27FOBFulQ7Zu4fdvdV7j5S//564CF3Pw94qP59z9q4bjlDg61T+6HB7LNDI2mPaDviyAuWjigUqb7ZWFC9Arik/vV3gH8CvjwLj5ObwSkpjMZ4pqS8R4Zo2qOsBcvR1cOM7z3MlkdeYNJ9RiWh2sQk0luKztwd+ImZbTezDfWx97r7AYD65/ck3dHMNpjZuJmNHzpUzkJeWkzOitWRbfdlbviJHHadVhKalWLRJiaR3lI0uK919w8CHwOuM7MPzfSO7n6bu4+4+8iiRdkleO1y2vzk8sC0cYjNWCM596h3viP52tPGIZZi0SYmkd5SKLi7+4v1zweBe4ELgZfM7AyA+ueDRS+yXdI2YGZ1/F2QEvjTxqHcUsNXE2rcs8Yh9oJVxSMKRaosnHM3s3cCA+7++/rXfwH8d+A+4Frg5vrnH7XjQtshEgjfmDLDnW4cYvX0kH8nLNReZJKO1Mt68TkzZYfqdCkWbWIS6R1FZu7vBf7ZzHYCjwIPuPv/pRbUP2pmvwE+Wv++K0Rm4UePn8g1HhUtNYy8G1GKRaT6wjN3d38euCBh/GXgI0UuarZEZuFlibYfSGo9kDUO6hMj0g/6qrdMZBYe6d0yaJbYYTGrpW601DDyWKAUi0jV9VX7gYjNl69gaEp/gqEBy+zdcvG5yZty08YhXmpYxXNhRaQ4BfdpjK4e5parL2ipErnl6gsyZ72/PvD7XOMQL5+MHLwhItXXV2mZqLwpjKTqlaxxiJdPauYuIkkU3GcgUqKYV9k59zL+TSLSOQru04gc1jF30Hhz8tSAO3cw+4zSSO15ZOZe5qlPItIZfZVzTztGLut4uchW/aTAnjUOse6TUN6/SUR6S18F98iiZZndECenNCib+n2Sbv83iUhn9FVwjyxaltUN8cb7d5/Sg/2E18azdPO/SUQ6p6+Ce6TnS1lb9SMVNhCbhav9gEj19VVwj9SEj64e5qo1wyd/ZiYHW6Stm2asp4bNnzuYaxzU4VGkH/RVtUy0siTpYIuRsxemBsO0ddOM9dTwKVFH30zui5M23qD2AyLV1lcz94iyKksip0QVuZ+IVJuC+zTKqiyJthFQ+wERSaLgPo2yKkuibQTWX7Q417iI9AcF92mUVVkSOd8V4KbRlXz64iUtC76fvngJN42ubOv1iUhv6asF1eGULf5ZuznLOtgicqJSw02jKxXMRaRFX83cu7m+O3K+q4hImr6auUdm4WU12Zo/d5A/JJQvZtWri4ik6avgDvnru6Nnm+YVrVcXEUkSTsuY2WIz+7mZPW1mu83sb+rjm81swsyeqH98vH2XW75IKWTa+apZ566qXl1E2qlIzv0t4Evu/mfAxcB1ZnZ+/bZvuPuq+sePC19lB0VKISPnrqpeXUTaKRzc3f2Auz9e//r3wNNA5fazRxZhI+euql5dRNrJvA1nbZrZUuBh4P3A3wJ/DbwGjFOb3b+ScJ8NwAaAJUuWrNm7d2/h65gtm8Z2seWRF5h0Z9CM9RctnpXSw09961f88rnDJ79fu2wh3/0v/67tjyMi1WBm2919JOm2wqWQZvbHwFbgi+7+GvBNYBmwCjgA3Jp0P3e/zd1H3H1k0aL0gyU6La1x2NiOibY/zqN7Wl8DH93zyoweZ2zHBGtv/hnnXP8Aa2/+WduvTUR6T6HgbmZD1AL7d939HgB3f8ndJ939BPAt4MLil9k5ZTUOu/H+3Ryf0jby+KRPe1hHo1Rz4sgxnLdLNRXgRfpbkWoZA74NPO3uX28aP6Ppxz4BPBW/vM4rq3FY9LAOnYcqIkmK1LmvBT4D7DKzJ+pjfwesN7NV1Kr49gCfK3SFHXZmSsuCbjmSTuehikiScHB3938m+SyJni59nGrjuuUtO1RhdloWLJg3xJGEVgNZtfHQ/S8+ItIZfdVbJqKsI+kitfHQ3f1yRKRz+q79QEQZR9JFu0+W1bVSRHpLW+rcixoZGfHx8fFOX4aISE+Z1Tp3ERHpPgruIiIVpOAuIlJBCu4iIhWk4C4iUkFdUS1jZoeAIm0hTwf+tU2X08v0PNToeajR81BT5efhbHdP7LzYFcG9KDMbTysH6id6Hmr0PNToeajp1+dBaRkRkQpScBcRqaCqBPfbOn0BXULPQ42ehxo9DzV9+TxUIucuIiKtqjJzFxGRJgruIiIV1NPB3cwuNbNnzey3ZnZ9p6+nU8xsj5ntMrMnzKyv2mua2e1mdtDMnmoaW2hmPzWz39Q/n9bJayxDyvOw2cwm6n8XT5jZxzt5jWUws8Vm9nMze9rMdpvZ39TH++5vomeDu5kNAv8L+BhwPrXj/c7v7FV11IfdfVUf1vPeAVw6Zex64CF3Pw94qIZsZHcAAAH1SURBVP591d3Bqc8DwDfqfxer3L1Sp6SleAv4krv/GXAxcF09LvTd30TPBnfgQuC37v68u78JfB+4osPXJCVz94eBw1OGrwC+U//6O8BoqRfVASnPQ99x9wPu/nj9698DTwPD9OHfRC8H92Hghabv99fH+pEDPzGz7Wa2odMX0wXe6+4HoPY/O/CeDl9PJ33BzJ6sp20qn4poZmZLgdXAI/Th30QvB/ekw7n7ta5zrbt/kFqK6joz+1CnL0i6wjeBZcAq4ABwa2cvpzxm9sfAVuCL7v5ap6+nE3o5uO8HFjd9fxbwYoeupaPc/cX654PAvdRSVv3sJTM7A6D++WCHr6cj3P0ld5909xPAt+iTvwszG6IW2L/r7vfUh/vub6KXg/tjwHlmdo6ZzQWuAe7r8DWVzszeaWbvanwN/AXwVPa9Ku8+4Nr619cCP+rgtXRMI5jVfYI++LswMwO+DTzt7l9vuqnv/iZ6eodqvbTrfwKDwO3u/tUOX1LpzOxcarN1gDnA9/rpeTCzLcAl1Nq6vgR8BRgDfgAsAfYBV7t7pRcbU56HS6ilZBzYA3yukXeuKjP798AvgF3Aifrw31HLu/fX30QvB3cREUnWy2kZERFJoeAuIlJBCu4iIhWk4C4iUkEK7iIiFaTgLiJSQQruIiIV9P8Bx6OTEu49/yYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from matplotlib import pyplot as plt\n",
    "\n",
    "\n",
    "#画图\n",
    "def draw(pred=None):\n",
    "    plt.scatter(data[:, 1], data[:, -1])\n",
    "\n",
    "    if pred is not None:\n",
    "        plt.scatter(data[:, 1], pred)\n",
    "\n",
    "\n",
    "draw()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Node col=0,value=1.00\n",
      "Leaf value=1.00\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "1"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#创建节点和叶子对象,用来构建树\n",
    "class Node():\n",
    "    def __init__(self, col, value):\n",
    "        self.col = col\n",
    "        self.value = value\n",
    "        self.lt_node = None\n",
    "        self.gt_node = None\n",
    "\n",
    "    def __call__(self, data):\n",
    "        if data[self.col] > self.value:\n",
    "            return self.gt_node(data)\n",
    "        return self.lt_node(data)\n",
    "\n",
    "    def __str__(self):\n",
    "        return 'Node col=%d,value=%.2f' % (self.col, self.value)\n",
    "\n",
    "\n",
    "class Leaf():\n",
    "    def __init__(self, value):\n",
    "        self.value = value\n",
    "\n",
    "    def __call__(self, data):\n",
    "        return self.value\n",
    "\n",
    "    def __str__(self):\n",
    "        return 'Leaf value=%.2f' % self.value\n",
    "\n",
    "\n",
    "node = Node(0, 1)\n",
    "leaf = Leaf(1)\n",
    "\n",
    "print(node)\n",
    "print(leaf)\n",
    "\n",
    "node.lt_node = leaf\n",
    "\n",
    "node(np.array([0, 0]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "836434.459201954"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#衡量一个数据子集y的一致性\n",
    "def get_error(data):\n",
    "    return np.var(data[:, -1]) * len(data)\n",
    "\n",
    "\n",
    "get_error(data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "((387, 3), (13, 3))"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#切分数据\n",
    "def split_data(data, col, value):\n",
    "    data_gt = data[data[:, col] > value]\n",
    "    data_lt = data[data[:, col] <= value]\n",
    "\n",
    "    return data_gt, data_lt\n",
    "\n",
    "\n",
    "data_gt, data_lt = split_data(data, 1, 0.2)\n",
    "data_gt.shape, data_lt.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<__main__.Node at 0x7fe9258b66d8>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def get_node(data):\n",
    "    min_error = np.inf\n",
    "    min_col = 0\n",
    "    min_value = 0\n",
    "\n",
    "    #遍历所有列\n",
    "    for col in range(M - 1):\n",
    "\n",
    "        #遍历所有值\n",
    "        for value in np.unique(data[:, col]):\n",
    "\n",
    "            #切分数据\n",
    "            data_gt, data_lt = split_data(data, col, value)\n",
    "\n",
    "            #如果一个子集是空,则说明另一个一定是全集,这种切分是没有意义的\n",
    "            if len(data_gt) == 0 or len(data_lt) == 0:\n",
    "                continue\n",
    "\n",
    "            #计算切分后的loss是否下降了\n",
    "            error = get_error(data_gt) + get_error(data_lt)\n",
    "            if error < min_error:\n",
    "                min_col = col\n",
    "                min_value = value\n",
    "                min_error = error\n",
    "\n",
    "    #前剪枝\n",
    "    #切分之后loss下降要超过1才切分,这是为了避免过拟合,否则定义为叶子节点\n",
    "    #这里故意训练一个复杂的树,方便后面后剪枝\n",
    "    if get_error(data) - min_error < 1:\n",
    "        return Leaf(data[:, -1].mean())\n",
    "\n",
    "    #创建新的节点\n",
    "    node = Node(min_col, min_value)\n",
    "\n",
    "    #切分数据\n",
    "    data_gt, data_lt = split_data(data, node.col, node.value)\n",
    "\n",
    "    #前剪枝\n",
    "    #切分后的子集不能太小,这是为了避免过拟合\n",
    "    if len(data_gt) < 3 or len(data_lt) < 3:\n",
    "        return Leaf(data[:, -1].mean())\n",
    "\n",
    "    #左右节点继续递归\n",
    "    node.gt_node = get_node(data_gt)\n",
    "    node.lt_node = get_node(data_lt)\n",
    "\n",
    "    return node\n",
    "\n",
    "\n",
    "root = get_node(data)\n",
    "root"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "---- Node col=1,value=10.00\n",
      "-------- Node col=1,value=17.00\n",
      "------------ Node col=1,value=19.00\n",
      "---------------- Node col=1,value=21.00\n",
      "-------------------- Node col=1,value=22.00\n",
      "------------------------ Leaf value=171.88\n",
      "------------------------ Leaf value=168.67\n",
      "-------------------- Node col=1,value=20.00\n",
      "------------------------ Leaf value=162.10\n",
      "------------------------ Leaf value=163.85\n",
      "---------------- Node col=1,value=18.00\n",
      "-------------------- Leaf value=154.02\n",
      "-------------------- Leaf value=151.35\n",
      "------------ Node col=1,value=13.00\n",
      "---------------- Node col=1,value=15.00\n",
      "-------------------- Node col=1,value=16.00\n",
      "------------------------ Leaf value=144.60\n",
      "------------------------ Leaf value=141.99\n",
      "-------------------- Node col=1,value=14.00\n",
      "------------------------ Leaf value=136.11\n",
      "------------------------ Leaf value=135.21\n",
      "---------------- Node col=1,value=12.00\n",
      "-------------------- Leaf value=127.37\n",
      "-------------------- Node col=1,value=11.00\n",
      "------------------------ Leaf value=117.77\n",
      "------------------------ Leaf value=113.34\n",
      "-------- Node col=1,value=6.00\n",
      "------------ Node col=1,value=7.00\n",
      "---------------- Node col=1,value=8.00\n",
      "-------------------- Node col=1,value=9.00\n",
      "------------------------ Leaf value=99.94\n",
      "------------------------ Leaf value=97.89\n",
      "-------------------- Leaf value=87.72\n",
      "---------------- Leaf value=73.28\n",
      "------------ Node col=1,value=0.00\n",
      "---------------- Node col=1,value=5.00\n",
      "-------------------- Leaf value=63.87\n",
      "-------------------- Node col=1,value=4.00\n",
      "------------------------ Leaf value=48.10\n",
      "------------------------ Node col=1,value=1.00\n",
      "---------------------------- Node col=1,value=3.00\n",
      "-------------------------------- Leaf value=34.23\n",
      "-------------------------------- Node col=1,value=2.00\n",
      "------------------------------------ Leaf value=38.08\n",
      "------------------------------------ Leaf value=39.65\n",
      "---------------------------- Leaf value=46.75\n",
      "---------------- Leaf value=79.76\n"
     ]
    }
   ],
   "source": [
    "#打印树的方法\n",
    "def print_tree(node, prefix=''):\n",
    "    prefix += '-' * 4\n",
    "    print(prefix, node)\n",
    "    if isinstance(node, Leaf):\n",
    "        return\n",
    "    print_tree(node.gt_node, prefix)\n",
    "    print_tree(node.lt_node, prefix)\n",
    "\n",
    "\n",
    "print_tree(root)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD4CAYAAAAXUaZHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3df5BU5Zkv8O8zTY9ANuXAMsnCBMSZZXFFAiNTyi00VyubxaQsHPGyC8Rc925qSapM3U1taiq613sFi5TW5ZrsH7ubLSwtc8uAisioZTYmZbLrjwqagQGBKFdBQAYWiAimwq9x5rl/dDf0j/OenvOc06f7nP5+qibMvDMnfWjGp99+3ud9XlFVEBFRurTU+waIiCh6DO5ERCnE4E5ElEIM7kREKcTgTkSUQuPqfQMAMGXKFJ05c2a9b4OIKFG2bdv2W1Vt9/peQwT3mTNnYmBgoN63QUSUKCJy0PU9pmWIiFKIwZ2IKIUY3ImIUojBnYgohRjciYhSqCGqZYiImk3/4BDWvbQXR06dxbS2CehbPBu93R2R/f8zuBMRxax/cAj3PrsLZ4dHAABDp87i3md3AUBkAZ5pGSKimK17ae/FwF5wdngE617aG9ljMLgTEcVs6NTZQOMWDO5ERDHLiAQat2DOnYgopKCLoyOOE/Bc4xYM7kSUSrWuRil+nKCLox1tEzxTMB1tEyK7L6ZliCh1CgF36NRZKC4F3P7Bocgfy7I42rd4NrKZ0hRMNiPoWzw7svticCei1ImjGqXAvDhanoGJLiMDgMGdiFLoiCOwusbDsCyOrntpL4ZHS6P58KiyFJKIyM80R+7aNR6GZXGUpZBERAZ9i2djQjZTMjYhm4k0p11gmbm3OL7lGrdgtQwRpU6hSiWOahnLzH3U8S3XuAWDOxGlUm93R02CebmMiGcgj3JDkgXTMkREIVhm7q6wH+XLAWfuRJRKcW1iahHvdIpf/twV9qOshuTMnYhSp39wCH3P7CzZxNT3zM6abGKy5M9dO1Fj3aEqIo+JyHER2V009pSI7Mh/HBCRHfnxmSJytuh7/xLZnRIRjdGaF/ZgeKSsjnxEseaFPXW6o1I3X9UeaNxiLGmZxwH8I4D/WxhQ1b8sfC4iDwM4XfTz+1R1flQ3SEQU1EdnhgONh9E2IYtTZyv/f9smZJ3X/PKdE4HGLarO3FX1FQAnvb4nIgLgLwBsjOyOiIgSZPWSORX59RbJjbvEsYM27ILqjQCOqeq7RWNXisgggI8B3Keqr3pdKCKrAKwCgBkzZoS8DSKi8KyLsBkRjBZVx1Qrg5zWNgGvnrkdxT+mCtw4cYv53suFDe4rUDprPwpghqp+KCILAPSLyBxV/bj8QlVdD2A9APT09ETcMoeI0iSOyhfruaZ+fWJc1716LhfYy18DXj13O0qz3HbmahkRGQdgKYCnCmOqel5VP8x/vg3APgB/EvYmiah5xdW+19pJ0pJiEVQGdpFo69zDlEL+GYB3VPVwYUBE2kUkk/+8E8AsAPvD3SIRNbO42vda8+CmJmUxFLpXTcuIyEYANwGYIiKHAdyvqo8CWI7KhdQvAHhARD4BMALgm6rquRhLRDQWcbXvneY4HalaJ8mbr2rHA9tvqMif/6+rXnNfFMMW1arBXVVXOMb/ymNsM4DN4W+LiCjHGnSDuvmqdjyx9ZDnuJ8HdtzgmT9/YMcNQK93/lzz/1P+gqARpma4Q5WIYtU/OIRFD/0CV97zIhY99IuqufM4NvwAwItvHQ00XiDqyJ/7pFi6zm3IBfOyj65zG4LethN7yxBRbAptAQq7RwttAQB3RYplw4+lU6N545MhxdI2MYvOM5WBfNJE98anoDhzJ6LYWNoCWHLuCzsnBRqPm6thpE8jycAY3IkoNpbZcZtjNusaB4ADH3oHftd4GALvs679cuenPdoV+I1bMC1DRA3tfFkZZLVxIDer39+6smLBsutUdDntgvvmv4YHtt9QEs1Vgf957WtY67gmjkViBnciMotj5+iZ4dFA4wCwb/xKz4XOfeNXIqodoAUb3/gAT2jli0bmjQ+wtneu5zXWypwgGNyJyMSyXd/SQdHCVcHit0moRYD3spWz/T8e9p/tW05i8qvMcb0gBMWcOxGZWHaOrl4ypyIXLfDvoBjHkXTApcBe/vFedmXEjxRPS2IGdyIysVSxDBw86bn4OHDQvZE9jiPpAHhuRPIaSwoGdyIysfRU2fjGB4HGAXdtul/NemFTULWxEnG9RYA7DRVleorBnYgC7xoFgL7FszEhmykZm5DNoG/xbOc1lvy05ZrOC947QDsv+OTP43qLAODWeVMDjVtwQZWoyVn7mPd2d2Dg4ElsfOMDjKgiI4I7FnT4XmPZOdrhKBv0O0x6kmEHaOG2Kvq9qP/kPbHH7BFRullb6vYPDuGpX39wMViPqOKpX3/gO+tfcf30QONA7h2C1zF2fu8QLLXxptk+0nvMHhElnDXQ+LUScM3eC2V+xbP9FddP9y3/2zRwCGUHHWFUc+Oux7HUxrdNyKLzbGUgH0se3HLMXq03MXHmTtTkLNv7AXs5X88Vk/FHl4+HAPijy8ej54rJvj//+j7vShrXuNXqJXOQLZuCZ1vEdwYO+B+z52JZrwiKM3eilAm6azSOJlbF9xa0K6SFa7+S33y68PhBd9x6zcD9xguPFXS9IigGd6IUsSyOWptYWRYSLakcAJ59Yvxy4dbCl97u4AHWskjcPziEzduGStYrNm8bQs8VkyML8EzLEKWIZXHUdAYobLtNLamcQmAv/9jf6t456qqk8auwsbKUasZxLiyDO1GKWBZHrflfy25TC2lx7Bz1iV5xnd4E2F5I4qiWYXAnShHLLLy3uwN3LOi4mEYYa/7XstsUyM2437/s0offDNzKemSeheXF0fpuKYiqwV1EHhOR4yKyu2hstYgMiciO/MdXir53r4i8JyJ7RWRxZHdKRFX1LZ6NbKas4iMjvoHGUq9e+Lkg44AtxWJJoMfRmKvA8uIYxzuLsczcHwdwi8f4D1R1fv7jJwAgIlcDWA5gTv6afxaRjMe1RFQrXrkSH5aj7wBbzxdLisXUJyZGrsVRvxfHhtihqqqvABhrEu02AE+q6nlVfR/AewCuC3F/RBSApebaOsu17Da16HLsHO2qsnM0LpbF0UbfofotEfmvAAYAfEdVPwLQAWBr0c8czo9VEJFVAFYBwIwZM0LcBhEVWGqurSy7TS2NFxXVWwCUi+tQEMAWqBt5h+oPAXQBmA/gKICH8+Ne/0aeb55Udb2q9qhqT3t79CvYRM2ovMdJtXEgXPvZoLtNR+GdYnE3BbClf6y7TS3dMS2Lo3HsUDUFd1U9pqojqjoK4BFcSr0cBlD8nuxzAI6Eu0Wi5nVf/y503fsTzLznRXTd+xPc17/L9+fLe7BUGwfCBcJ7n92FoVNnobi0YcovIN44fotniuXG8Vuc13S2Tww0DuQWOdctm4eOtgkQ5MoS1y2b57vIafn7ALZA3dvdgQeXzi25vweXzq3/DlURmaqqhZqi2wEUKmmeB7BBRL4PYBqAWQDeDH2XRE3ovv5dJYcoj6he/DqqczYB763wf3nd9KqBxi/X7Lr2yKmz6ERlikUuuFMY+0+cCTReEHS3qeXvU3icwvVB2hZYdsMGUTW4i8hGADcBmCIihwHcD+AmEZmPXMrlAIBvAICq7hGRpwH8BsAnAO5WVXePTSJy8qsjjzK4u0ohq22FP3LqrGdbgK5T7vy4JddsKbkEci+OQdYDwixy1jpQW4ylWmaFqk5V1ayqfk5VH1XVr6nqXFX9vKouKZrFQ1W/p6pdqjpbVf+1trdPlF6WoDYx6/2ftGscsJdC7hvvXbO+b7y7Zt1Sh2/JuRfe9RS/YD2x9ZBvWiuOjUVx4g5VopgEzZ9btI7z3lbiGgfspZACR82671UIXIdvKbm07J6NY5EzTuwKSRSDuPLn1g6PQbsu5n4IAerjcvzq8F1pjZ4rJmPD1kMlFTUt+XEXy7sea+68UXHmThQDy0zSko6wpBb2O9Ir+33SKwBMReuWvPa6l/ZWlEqO5sddLM9d2jC4E8XAMpO0XGPJaVvTKwpHWwCfaywvPpYXBEsqx1oK2agY3IliYJlJfqrVO0/uGi8YKUt7lH9dwXiyRdc5R1uAc+50jiWvbTkGcG3vXNy5cEZJM687F87wTYHF0WM9Tsy5E8VgxfXTS3LuxeMuv7/gXUXsGgdylS9eh0lXO+nIom1iFp1nKgP5JJ+g29vdgU0Dh0rOP712xuU1OQZwbe/cQOsZcfR7iRNn7kQxWNs7F4u6ShcAF3VNjnQxFbBVvli7LlqC7n39uyoOtn5930nfyiHrInFQLIUkosD6B4ew/dDpkrHth043RD6309F1sVq1jFdjLr9xwLawbEnLWLAUkogCs25tD6o1I3gns6KirPGqkY2+1wXtugjkmpF5pfP9mpRZFonPD3unoVzjVmkrhWRwJ4qBNZ8btP68ENjL12nfyawAcNrzGitLkzLLC8KZYe+eka7xMBqxjYAV0zJEMTDVnxuOpPMK7F5j9XLZOO+Q4xonOz6jRDGwnJnZ6IHa4pxjtu0aB8L1m29mDO5EMTCdmWk5tsjAupvTcnuWdzDWfvPNjsGdKCDLaT1xHX9nKWu0ttT96kLv4zFd44D9YIugB28QF1SJAilsUS9UvhS2qAPwDTYZEc9g6Tc71lEAZWmYQqB2XdV5YUNFTr6wCHsgwnsDbGeoNurBFmnE4E4UgLWk0TI7vnHiFrx65vaSMdXc+Os+9xi0rNE6cweC7wIFGKjjwuBOFIA1vdLhOIGoo8ohyrM2bSzpDZNpETwc8aYay71R42POnZpa0Py5dfHRkmseOHjSswnYwMGTjits0rYzk3I4c6emZcmfW1MYvd0dWPLc1ZAMLh5yoQBaut0bi+I6Q9XrgOw7FtQuddI/OJSaXaCNjDN3alqWFq+uVEXVFMbqy9GCS73TBfn/+FZf7rzE8kLi6sjo16mxf3AIm7cNlZw3unnbUE363qStZ3ojY3CnpmVpCWA5DMPKkgKydGqMs4952nqmN7KqwV1EHhOR4yKyu2hsnYi8IyJvicgWEWnLj88UkbMisiP/8S+1vHmiMMwtXgMe8GxlOU3I0qkxzj7maeuZ3sjGMnN/HMAtZWM/B3CNqn4ewP8DcG/R9/ap6vz8xzejuU2i6FlaAvgd8OzHctiR5TQhizj7mKetZ3ojq7qgqqqviMjMsrGfFX25FcB/ifa2iGrP0hLAWgo5CqBFKzckjQrgd2iepY48qL7Fs0sWloGxV8sEXRwN81gUTBQ5978G8K9FX18pIoMi8u8icqPrIhFZJSIDIjJw4oRPfw2iGrGkCFytaf1a1gLAH5/bgNGywzBGNTceJcv99XZ34MGlc0u29z+4dG7VChbL4qj1sSi4UKWQIvI/AHwC4Mf5oaMAZqjqhyKyAEC/iMxR1Y/Lr1XV9QDWA0BPT0+NspZEbm0Ts57Hz/md8GPpYQ7k0g5dpyoDedQbhS4b14KzHh0Wq7XUtewate7W5Q7VeJhn7iJyF4BbAXxVNbcWr6rnVfXD/OfbAOwD8CdR3ChR1KwHL1tYq2yCbrLyCux+42FwcbSxmWbuInILgO8C+M+qeqZovB3ASVUdEZFOALMA7I/kTomqCJr/jevg5QKv3aZ+LJusrE3ALKY52hZwcbQxjKUUciOAXwGYLSKHReTrAP4RwKcB/Lys5PELAN4SkZ0AngHwTVWNdq80kQdL/jeug5cBYM0LeypSN6OaG3ex1ISHaQIWFNsWNLaxVMus8Bh+1PGzmwFsDntTREFZ8r/WtEzQc00BeOb2/cYBW9ojziZgaTtQOm3YW4ZSwRIILWmZ4nNNy8ejPoDakvaIu9SQi6ONi+0HKBUsm2Ms11jPNZ2Y9f5PzTUO2E8tYqkhAZy5U0pYZqw3X9WOJ7Ye8hyPWuu4DM54VKy0jnNvYeKpRRQGgzulgqVt7Zbt3outW7YPuXeFGg+ttvR8ARioyY7BnRpS0LJGV9vanismO6/7/YWRQOOA7VxTIN4SRSKAwZ0akKW+27pbMmjli+UAaiDeEkUigMGdGpAlUFuqZayVL0EPoAZ4TinFj9Uy1HAsgTrOyhcLS3thojAY3KnhWAK1KXgaF0ctLO2FicJgcKeGY6nvfvGto4HGAVtsb5vgaFngGC9gky2KG4M7NRzLRhzL9v5RVLYaUM2Nu9w6b2qg8QKeQERx44IqNaQ46rsXtGzCtk+WlYypAgsymzDouMaaXuEJRBQ3BndKhYnZFs8doH7b+z86M4xOeFS+XIj+MGk22aK4MbhTQwq6icmyvd8iTA9zy7uRoM8DUQFz7tRw+geH0PfMzpLe7H3P7PTtzW7p8GhZHI2zh7mlRz1RAYM7NZw1L+zB8EjpSufwiPoebGFZsLQsjvZ2d+COBR0X2waMpYeNleWwDqICBneqqaBngAK2yhdLnbtlcdTVw6YWs2mWT1IYDO5UM3GmFSx17l65c79xIN7ZNMsnKQwGd6qZOAOhZbbv6sjo16kxztk0zyilMBjcqWYaPa1g6dQY52yapypRGFVLIUXkMQC3AjiuqtfkxyYDeArATAAHAPyFqn6U/969AL4OYATAf1fVl2py59TwwpQNBtU2Iet58IVf5YulxzrPKKWkGMvM/XEAt5SN3QPgZVWdBeDl/NcQkasBLAcwJ3/NP4tItIXGlBhxphVWL5mD/a0r8f5llz72t67E6iVznNdYZu6cTVNSVJ25q+orIjKzbPg2ADflP/8RgH8D8N38+JOqeh7A+yLyHoDrAPwqmtulJIlzV2bvc1dDW8qafkluHN3evdmtPdY5m6YksO5Q/ayqHgUAVT0qIp/Jj3cA2Fr0c4fzY9SkrIEw6AlJQGU3x2qde+M8IJsoblEvqHr99+T5HldEVonIgIgMnDjBntZ0SfEJScUf5cfbFXMlUvwOsWOPdUoza3A/JiJTASD/5/H8+GEA04t+7nMAjnj9H6jqelXtUdWe9nbOlOgS0wlJhuje6NU8RGFYg/vzAO7Kf34XgOeKxpeLyGUiciWAWQDeDHeLlGSWHaoWqt692f3On+YmIUqzsZRCbkRu8XSKiBwGcD+AhwA8LSJfB3AIwDIAUNU9IvI0gN8A+ATA3ao64vl/TIkTtENhYYdqoWywsEMVQOQLkrOGN+Ld7IqSMdXc+D7HNdayRnZqpCQYS7XMCse3vuj4+e8B+F6Ym6LGYwnUfjtU/YJhYbZdvqCq6l4kXXH9dHRurVxwvXPhdI+fRsl9N+oLFlEY7OdOY2IJ1Nac9jX6JHZjecmYam78N45req6YjA1bD5UckdeSH/cTtJrH+oJFFDcGdxoTS5OttolZz94ubRP9D5M+MzzqfUKSz+mm617aW/Hd0fx4lEGXi7CUFOwtQ2NiabL1+/OfBBoPw/LiY8FFWEoKBncaE8tW/Qsj3t9zjYfR4niNcY1bsVMjJQWDO42Ja0t+ta36cRl1vF64xq3YW4aSItE5d5akxSfuboiNjL1lKAkSG9wLhygXztosHKIMsCStFqxNwCw9YogovMQGd79DlBncayPojLW4R0z5OODdqRGw9WYnolKJzblbjlWjcIK2EjD1iEGuN3u2bCU02yK+vdmJqFRiZ+4Urzh3ZsbZB54orRjcaUzi3pnJRUuicBKblqF4WTYJWTo1Wrny8czTU7NKbHCfmPW+ddc4hWPZodp5YcPFYF78UYtqGebpiUoxLUNjYtmhmm3xDuS1eP1lnp6oVGKD+5lh7yZSrnEKx3KYtOufolb/RMzTE12S2OAOcINMnHiYNFGyJDZBvX+84xDl8e5DlMnuxbeOBhoH3AdrRNzLi4g8JHbmLvDeIEO1Ydk0ZjizmogiktiZOzU+S4UNEUUjucGd08JYWerILRU2RBSNxAb3ODfIpFHQPjG3zpsaaBzgzJ2onszBXURmi8iOoo+PReTbIrJaRIaKxr8S5Q0XxLlBJm0KfWKGTp2F4lKfGL8Ab1lQ5cydqH7MC6qquhfAfAAQkQyAIQBbAPw3AD9Q1f8TyR06ZEQ8AzlnhdVZ+sRYFlQnOQ7InlTlgGwiCi+qapkvAtinqgclpuDKWaGd9TDpoPsKXP8U/Cciqr2ocu7LAWws+vpbIvKWiDwmIpO8LhCRVSIyICIDJ06cCPyAcR2ITDnFB2+U7Ctode8rOO1x4IbfOBFFJ3RwF5FWAEsAbMoP/RBAF3Ipm6MAHva6TlXXq2qPqva0twff5RjXgciUYzl4Y5qjNYFrvFjQBV8iKhXFzP3LALar6jEAUNVjqjqiqqMAHgFwXQSPQQnUt3g2JmQzJWNjOVTbsuBLRKWiCO4rUJSSEZHi2rjbAeyO4DEqcGu7XcbxJLnGrXq7O/Dg0rnoaJsAQa7J2INL51Zt7uW34EtEYxNqQVVEJgL4EoBvFA3/bxGZj9x2ogNl34sM9zDZWVJahUXQ8gVVVf8XVEunxiOOhV3XOBFVCjVzV9UzqvqHqnq6aOxrqjpXVT+vqktU1V0IHYKrnI5ldtVZXhivzWzy3FdwbWaTz1U2YXL1RJST2MZhLLOL1+mzw+hEZdmjIPrKl77Fs0sO4wbGlqsnoksSG9ytZXb9g0M8rcdgmuOwjlrMpnmqElF4iQ3ulmBTqMIozAgLVRgAGDiq6Fs8G9/ZtBMjRYn5TIvUbDbNU5WIwkls4zBLmR2rMHJcZ5j6nW06cPBkSWAHgJFRxcDBkxHeGRFFJbEzd8tbd1Zh5GQzLdg7bnlF5cs1+qTzmo1vfOAcX9s7N+pbJKKQEhvcgeBv3ePMG8cp6DrCblnuubt0N5YDOO15DXv5ECVLYtMyFtYdk42sf3AIfc/sLNnN2ffMTt/dnJZWAkSULE0V3K07JhvZmhf2YHikdPY8PKJY88KeOt0RETWCpgruQG5h8D9On4MC+I/T5xK/IGjps05E6ZfonHtQ9/XvwhNbD138ekT14tfNtChobSVARMnRVDN3v4qPpLIcXG05otDyOERUP001c09jxcecaZ/G6/sqU0tzpn3a97qgZ82uXjIHf/f0jpLmYi2SGyeixtNUM/c02rr/o0DjALCoa3Kg8YLy82l5Xi1R42JwTzjLu5EDH3pv2nKNA7nNYsNlO1SHR7XpdvcSJUVTBfcOx2Yl13haWQ7I5u5eomRpquCexk1MFq50il+ahT3WiZKlqYJ7b3cH7ljQcTGIZURwx4LadB+0HPAc16HQllQOXxiJkqWpqmX6B4ewedvQxSA2oorN24bQc8Vk3wAftHeLpbVwnO2IOxw9dvzSU+yxTpQsTRXc/Vr+Rhl0LY9juQbI1Zmf8jigxK/+/Oar2ks2cxWP+2GPdaLkaKq0jGUh0dID3rL4aF2wXL1kTsWuUoF//fkv3zkRaJyIkqepgnuLY73QNQ7Ygq5l8dG6YDlw8GTFwdaaH3exvpDEtSZAROGFCu4ickBEdonIDhEZyI9NFpGfi8i7+T8nRXOr4Y061gtd44At6LrSG35pD+uC5YY3KtMrfuOA7e9USE8Vtxa+99ldDPBEDSqKmfvNqjpfVXvyX98D4GVVnQXg5fzXidW3eDaymdKpfTbjf3aoJe1hbUdsecHiEYVE6VeLBdXbANyU//xHAP4NwHdr8DiBCVCRwiiM+/LKe/iwpj2sC5b7W1dWdHj06x3T292BgYMnsfGNDzCiOqaSUG5iIkqWsDN3BfAzEdkmIqvyY59V1aMAkP/zM14XisgqERkQkYETJ+JZyHPFZL9Ybdl2H+eGn0JgL//Y37rSeY2rJNQvxcJNTETJEja4L1LVawF8GcDdIvKFsV6oqutVtUdVe9rb/UvwojJpYhb7W1fi/csufexvXYlJE91lg5YZqyXnbiUtjiPzfP5lLSkWbmIiSpZQwV1Vj+T/PA5gC4DrABwTkakAkP/zeNibjMq20WWes9xto8uc17Q5Ar9rHIi51NDwdsTygpXGIwqJ0syccxeRTwFoUdXf5T//cwAPAHgewF0AHsr/+VwUNxoFUe9Zrl8gPF82w602Dtjq6YHgO2EBuBcMfBYSpjl2qFZLsXATE1FyhJm5fxbAayKyE8CbAF5U1Z8iF9S/JCLvAvhS/uvGYAiEZ4ZHA41bWUsNFZeOzbs4pv7rCEyxEKWfeeauqvsBzPMY/xDAF8PcVM0ovAN5AxzEZG0/0HluQ8XiqSrQeX4DDjiuYZ8YovRrqt4yloOhLb1bMiKeHRb9WupaSw0zIp5lj9VOSWKKhSjdmqr9gOVg6NVL5iBb1p8g2yK+vVsWdnpvynWNA/ZSwzSeC0tE4TXVzB0IfjC0JYXxm6O/CzQO2Ds1Wt4lEFH6NV1wt+h9fh56dRgYD+AcgOezQPdvnT//0ZnKNI7fOGAvn+TMnYi8MLhXs2YKVIdLcvKqw5A1U4D73QE+qDA5d8vM3VR2SUSJweBeRXlgB/I9ajzGC1ozgncyKyoWbq8a2eh8HGvtuWXmHuepT0RUH021oOo6Rs7veDmLQmAv/3gns8J5jaX7JGD7O7HDI1H6NVVwN/V8MWzvLwTzamPlRsoalJV/7cXyd2KHR6L0a6rgblm0PC8Zzx2g5yXjfYHRmhf2VPRgH9XcuB/L34kdHonSr6mCu6Xny09v24Xz2lJSF39eW/DT23ZFem8fnRn27FjpV2ED2GbhbD9AlH5NtaBqqSzp7e7AfQdfKTnYYsX107HWZ+HRshO2uC97+Thw2vlYE1sz+P2FyiZmE1vd7yzYfoAo/ZoquFsrS7wOtui5YrIzGHZecPR7ueDu9+LK01dzxiOw+40XsP0AUbo1VXC3MDf0CrgT1spyuhQRpV9T5dwtGr2yxJVSYvsBoubG4F5FXJUlhZx8tbFyK66fHmiciJoDg3sVcVWWLBi3ybNj5YJxm3yvW9s7F3cunHFxpp4RwZ0LZ2Bt79xI74+IkqWpcu4dji3+frs546oscbUebhtDOf3a3rkM5kRUoqlm7o1c333a40AQv3EiIj9NNXO3zMLjarJlqVcnInJpquAOBK/vtpZCBmWtVyci8mJOy4jIdBH5pYi8LSJ7RORv8+OrRWRIRHbkP74S3Sx/qjkAAAUvSURBVO3Gz1IK6Tpf1e/cVdarE1GUwuTcPwHwHVX9UwALAdwtIlfnv/cDVZ2f//hJ6LusI0sppOXcVdarE1GUzMFdVY+q6vb8578D8DaA1O1ntyzC9nZ3YN2yeehomwBBrhpn3bJ5vmkc1qsTUZREIzhrU0RmAngFwDUA/g7AXwH4GMAAcrP7jzyuWQVgFQDMmDFjwcGDB0PfR63c17+rsnFYDUoPv/rIr/D6vpMXv17UNRk//pv/FPnjEFE6iMg2Ve3x+l7oUkgR+QMAmwF8W1U/BvBDAF0A5gM4CuBhr+tUdb2q9qhqT3u7z2EZdeZqHNY/OBT547x5oPQ18M0DH43pcfoHh7DooV/gyntexKKHfhH5vRFR8oQK7iKSRS6w/1hVnwUAVT2mqiOqOgrgEQDXhb/N+onrSLo1L+zB8Ejpu6jhEa16WEehVHPo1FkoLpVqMsATNbcw1TIC4FEAb6vq94vGpxb92O0Adttvr/7iahzmOpSj2mEdPA+ViLyEqXNfBOBrAHaJyI782N8DWCEi85Gr4jsA4Buh7rDOpjlaFjTKkXSN3rWSiOrDHNxV9TV4HyyU6NLHcn2LZ5fsUAVq07KgbUIWpzxaDfjVxgON/+JDRPXRVL1lLHq7O/Dg0rklZY0PLp0beeMwS2080Nj9coiofpqu/YBFHEfSWbtP8jxUIvISSZ17WD09PTowMFDv2yAiSpSa1rkTEVHjYXAnIkohBnciohRicCciSiEGdyKiFGqIahkROQEgTFvIKQB+G9HtJBmfhxw+Dzl8HnLS/DxcoaqenRcbIriHJSIDrnKgZsLnIYfPQw6fh5xmfR6YliEiSiEGdyKiFEpLcF9f7xtoEHwecvg85PB5yGnK5yEVOXciIiqVlpk7EREVYXAnIkqhRAd3EblFRPaKyHsick+976deROSAiOwSkR0i0lTtNUXkMRE5LiK7i8Ymi8jPReTd/J+T6nmPcXA8D6tFZCj/e7FDRL5Sz3uMg4hMF5FfisjbIrJHRP42P950vxOJDe4ikgHwTwC+DOBq5I73u7q+d1VXN6vq/Cas530cwC1lY/cAeFlVZwF4Of912j2OyucBAH6Q/72Yr6qpOiXN4RMA31HVPwWwEMDd+bjQdL8TiQ3uAK4D8J6q7lfVCwCeBHBbne+JYqaqrwA4WTZ8G4Af5T//EYDeWG+qDhzPQ9NR1aOquj3/+e8AvA2gA034O5Hk4N4B4IOirw/nx5qRAviZiGwTkVX1vpkG8FlVPQrk/mMH8Jk63089fUtE3sqnbVKfiigmIjMBdAN4A034O5Hk4O51OHez1nUuUtVrkUtR3S0iX6j3DVFD+CGALgDzARwF8HB9byc+IvIHADYD+Laqflzv+6mHJAf3wwCmF339OQBH6nQvdaWqR/J/HgewBbmUVTM7JiJTASD/5/E6309dqOoxVR1R1VEAj6BJfi9EJItcYP+xqj6bH26634kkB/dfA5glIleKSCuA5QCer/M9xU5EPiUiny58DuDPAez2vyr1ngdwV/7zuwA8V8d7qZtCMMu7HU3weyEiAuBRAG+r6veLvtV0vxOJ3qGaL+36BwAZAI+p6vfqfEuxE5FO5GbrADAOwIZmeh5EZCOAm5Br63oMwP0A+gE8DWAGgEMAlqlqqhcbHc/DTcilZBTAAQDfKOSd00pEbgDwKoBdAEbzw3+PXN69uX4nkhzciYjIW5LTMkRE5MDgTkSUQgzuREQpxOBORJRCDO5ERCnE4E5ElEIM7kREKfT/AQiW1r+anjd+AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#预测\n",
    "pred = np.empty(N)\n",
    "for i in range(N):\n",
    "    pred[i] = root(data[i])\n",
    "\n",
    "draw(pred)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
